Generalization of the singlet sector valence bond loop algorithm to antiferromagnetic 

ground states with total spin Stot = 1/2 
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We develop a generalization of the singlet sector valence bond basis projection algorithm of Sand- 
vik, Beach, and Evertz (A. W. Sandvik, Phys. Rev. Lett. 95, 207203 (2005); K. S. D. Beach and 
A.W. Sandvik, Nucl. Phys. B750, 142 (2006); A. W. Sandvik and H. G. Evertz, arXiv:0807.0682, 
unpublished.) to cases in which the ground state of an antiferromagnetic Hamiltonian has total 
spin Stot ~ 1/2 in a finite size system. We explain how various ground state expectation values 
may be calculated by generalizations of the estimators developed in the singlet case, and illustrate 
the power of the method by calculating the ground state spin texture and bond energies in a L x L 
Heisenberg antiferromagnet with L odd and free boundaries. 
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INTRODUCTION 

Understanding the ground states of strongly corre- 
lated condensed matter systems is a central problem 
in computational physics. Several approaches have had 
a degree of success in this endevour. These include 
various sophisticated quantum monte-carlo techniques 
for sampling the partition function Tr (exp(— TJ/T)) of 
a system with Hamiltonian H at temperature T, and 
using this sampling procedure to estimate the ther- 
mal expectation values of various operators (0)t — 

Tr (dexp{^H/T)j /Tr{eiip{-H/T)) Although 

these can be used for relatively large finite-size systems, 
they are intrinsically finite-temperature methods and ac- 
cessing the very low temperature regime involves doing 
calculations at successively lower temperatures and then 
extrapolating. 

Other approaches include various exact diagonaliza- 
tion techniques that obtain the lowest energy state in a 
given sector. These are severally constrained by memory 
requirements in terms of the system sizes they can han- 
dle. While this problem can be overcome in one dimen- 
sion by the sophisticated density matrix renormalization 
group method Q , there is as yet no generalization of this 
method that works equally well in higher dimension, al- 
though there has been considerable progress recently [1] . 

Recently, Sandvik and collaborators have developed 
an extremely elegant and sophisticated projection algo- 
rithm [iHl] that essentially solves the problem of cal- 
culating the ground state expectation values of quan- 
tities in a large class of S" = 1/2 antiferromagnetic 
spin systems which have ground state in the total spin 
Stot — singlet sector. This singlet sector algorithm 
works in the 'bipartite valence-bond' basis for singlet 
states (see below), and exploits the over-completeness 
of this basis to develop a procedure Q for evaluat- 
ing the expectation values of some observables O in 
the (unnormalized) ground state Itpg) = (— i?)™|s) ob- 
tained by acting on an arbitrary singlet state \s) with 



a large power of —H. The key to its success is an 
extremely efficient Q procedure for stochastically sam- 
pling (s|(-iJ)™0(-i?)™|s)/(s|(-i7)2™|s), which allows 
one to handle large systems with as many as 10'* spins in 
favorable cases. 

The ground state spin of a finite system made up of 
spin-half variables interacting antiferromagnetically nat- 
urally depends on the nature of the finite sample: if the 
total number of spin-half variables is even, one expects a 
ground state in the singlet sector, while systems with an 
odd number of spins will have a ground state in the total 
spin Stot = 1/2 sector. For instance, an Lx L square lat- 
tice Heisenberg antiferromagnet with periodic boundary 
conditions and L even will have a singlet ground state, 
while the same magnet with L odd and free boundaries 
will have a ground state spin of Stot = 1/2. In many situ- 
ations, it is useful to be able to handle both kinds of finite 
systems. For instance, if one wants to model experiments 
that dope insulating antiferromagnets with non-magnetic 
ions like Zn [13] that substitute for the magnetic mo- 
ments, it is convenient to study L x L periodic systems 
with L even as before, but with one spin removed from 
the system to model the missing-spin defect introduced 
by Zn doping. 

The original valence-bond projector loop algorithm P, 
Q allows one to study the singlet sector ground states 
of systems with an even number of spins. Here we ask 
if it is possible to come up with an analogous procedure 
in the total spin Stot = 1/2 sector of systems with an 
odd number of spins in order to compute properties of 
the Stot = 1/2 doublet ground state of antiferromagnetic 
systems with an odd number of spin-half variables inter- 
acting antiferromagnetically. As our results demonstrate, 
the answer turns out to be very satisfying: Using a ju- 
diciously chosen basis for the Stot = 1/2 sector of such 
systems, we find that is indeed possible to construct an 
analogous procedure that works as well in the total spin 
Stot = 1/2 sector as the original singlet sector algorithm 
of Sandvik and collaborators [ll-|3| • Here we detail several 
aspects of this generalization. To illustrate the power of 
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the method, we also show results for the ground state 
'spin texture' in the S^^^ = Stot = 1/2 ground state of an 
Lx L square lattice S* = 1/2 Heisenberg antiferromagnet 
with open boundary conditions and L odd and as large 
as L = 97. 
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BASIS 

A judicious choice of basis is the key to generalizing 
the original singlet sector valence-bond projector loop 
QMC algorithm to the study of bipartite spin-half an- 
tiferromagnets with Nb B-sublattice sites, Na — Nb + 1 
A-sublattice sites, and a doublet ground state in the 
Stot = 1/2 sector. While other choices may also be pos- 
sible, we find it convenient to use the basis 
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FIG. 1: Action of singlet Projector (P, denoted by dashed 
line) in valence bond basis. 



{|Aa/a)} = {\[A{b,)b,],[A{b2)b2]...[A{bMs)];afa}} 

(1) 

Each member of this basis has one A-sublattice spin S^^. 
in either the |a/cr =t) = l^^^. = +1/2) or the |a/(T =|) = 
\S^^ — —1/2) state along the quantization axis z, while 
the Nb spins Sfc^ on the B-sublattice sites each form a 
singlet state ('valence-bond') 

Uib^m . (i-^(Mt,^a)H^(M;,^.t)) 

V 2 

with a partner S_4(h.) on the A-sublattice. All basis states 
are obtained by allowing all possible a/, two choices for 
a, and all possible 'matching' functions A consistent with 
a given choice of 'free spin' af. Note that this basis set 
is actually a union of two distinct basis sets 

{|Afl/t>} = {\[Aib,)b,],[A{b2)b,]...[AibN,)];aft)} 

(3) 

and 

{\A,afi)} = {\[A{b,)b,],[A{b2)b,]...[A{br,,)];af i)} 

(4) 

corresponding to the two allowed choices for the con- 
served quantum number S^^^ in the Stot = 1/2 sector of 
an SU{2) invariant composed of an odd number of spin- 
l/2s. 

This basis is (over-) complete in a manner entirely anal- 
ogous to the bipartite valence-bond basis that was used 
in the original singlet sector algorithm 0,11]. This may be 
seen as follows: Consider adding one extra i3-sublattice 
site b^B+i to our system to make the total number of 
spins even. The singlet sector of this larger system is 
spanned by the (over-)complete bipartite valence bond 
basis. States in this basis are in one-to-one correspon- 
dence with possible pair-wise matchings V that 'find' a 



A-sublattice 'partner' ■p(&i) for each B-sublattice site bi 
to form a singlet: 

\V) EE \[Pibl)b^],[Pib2)b2]...[PibN, + lbN, + l)]) (5) 

Now, by the laws of angular momentum addition, sin- 
glet states of the larger system can only arise from ten- 
sor products of the additional spin-half variable at site 
bNt+i with the Stot = 1/2 states of the smaller system. 
Therefore, to check for (over-)completeness of our pro- 
posed basis for the smaller system, we only need to check 
whether all states in the bipartite valence bond basis of 
the larger system are obtainable as tensor products of 
states of the additional spin S^^^^^ with states in our 
proposed Stot = 1/2 basis. This is certainly the case, as 
is readily seen by identifying af with V{biq^j^i) and A 
with the restriction of V to the domain (6i, 62, ■ • ■ bN,,)- 
Our proposed basis is thus overcomplete in a manner en- 
tirely analogous to the original bipartite valence bond 
basis for the singlet sector. 

In practice, for SU(2) symmetric spin Hamiltonians of 
interest to us here, we will additionally exploit the con- 
servation of the z component of spin and restrict atten- 
tion to the basis {|[A(5i)5i], [^(62)62] . . . [A{bN^)];af t)} 
that only spans the S^^f. — Stot — 1/2 sector a system 
with Na — Nb + 1 spin-l/2s on the A-sublattice, and 
Nb spin-l/2s on the B sublattice. 

OVERLAPS AND OPERATORS 

We now indicate the changes that arise in the formu- 
lae for the wavefunction overlaps between basis states, 
and for the action of exchange operators when working 
in the Stot = 1/2 sector. As is well-known, the wavefunc- 
tion overlap between two bipartite VB basis states \V) 
and IV') of the singlet-sector basis are determined by the 
number of cycles needed to go from the permutation V 
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to the permutation V' . More pictorially, one may con- 
sider the overlap diagram of the two valence-bond covers 
viewed as 'complete dimer covers' or 'perfect' matchings. 
This overlap diagram contains A^ioop (closed) loops of var- 
ious lengths l^, such that each site is part of exactly one 
loop (see Fig[T|). Knowing this overlap diagram, one may 
calculate the corresponding wavefunction overlap to be 
{V\r') = 2^i°°p/2^=/2, where Ns, the total number of 
spins is assumed even. 

Generalizing this to states in our basis for the Stot — 
1/2 sector, we note that the corresponding picture is now 
terms of the overlap diagram of two partial valence bond 
covers {A, a/) and {A! , a'j), each of which leaves one site 
free (uncovered by a valence bond). Such an overlap 
diagram necessarily involves exactly one 'open string' 
of length If connecting a/ to a'p in addition to iVioop 
(closed) loops of various lengths la (see Fig [J). An ele- 
mentary calculation reveals that the wavefunction over- 



lap in our Stot = 1/2 case is given as {Aafa\A' a'jcr') = 
(5^^'2^'°°p/2(^=~i)/2, with Ns, the number of sites, now 
taken to be odd. 

The original singlet sector algorithm relies heavily [iHl] 
on a particularly simple action of operators Pij ~ rjirjjSi ■ 
Sj -1-1/4 on basis states {V). Here r/i = +1 {jji — — 1) 
for i belonging to the A-sublattice (B-sublattice), and 
thus, the operator Pa^b^ that connects an yl-sublattice 
site Ua to a i?-sublattice site bp is precisely the projection 
operator that projects to the singlet state of the two spins 
Sa„ and Sb/j. Our key observation, which allows us to 
generalize this algorithm to the the Stot = 1/2 case, is 
that the action of Pij on states in our 5*40* = 1/2 basis 
remains simple. This is seen as follows: If neither i nor 
j correspond to the 'free' spin, Pij acts exactly as in the 
earlier singlet sector case (Fig[Ha,b)): 
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I 

On the other hand, if either i or j correspond to the (Fig [He, d)); 
free site a / , one can easily check that the following holds 



PafbJ-[A{bo,)bo,]... ■,af t) 
PafAibJ\-[Aiba)ba]... ■,af t) 

Thus Pij either causes no change or rearranges exactly 
one pair of valence bonds to give a new basis state with 
amplitude 1/2, or reconnects one valence bond to move 
the free spin to give a new basis state, again with ampli- 
tude 1/2. The important thing to note is that these rules 
are in complete analogy to the original singlet sector case. 

By analogy to the original singlet sector work [iHl], 
this allows us to formulate a convenient prescription 
for the calculation of {A'a'f t |Pij|^a/ t) between two 
of our basis states by writing {A'a'f f \Pij\^o,f f) = 



= -|...[a/,6J... ;A{ba)t), 
= i|...[a/,6J... ;^(&a)t)- 

(7) 

I 

Wij {A'a'f t \Aaf t) and developing rules for the weight 
Wij by comparing the overlap diagram of (yl'a'j and 
{Pij \Aaf t) with the original overlap diagram of {A'a'f f) 
and {Aaf f): If the action of Pij makes no changes 
in the original overlap diagram, Wij — 1. In addition, 
Wij = 2x1/2 = 1 if a loop is split into two loops or 
the open string is split into one loop and another open 
string; here, the factor of two comes from the fact that the 
number of loops in the overlap diagram increases by one, 
while the factor of half has its origins in the reconnection 
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FIG. 2: Action of a singlet projector on the overlap loops 
for closed loops (a, b) and open string (c,d). Corresponding 
weights are given below within parenthesis. 



amplitude of one-half in Eqn[7l On the other hand, if two 
loops fuse into one, or if the open string fuses with a loop 
to give a larger open string, then Wij = (1/2) x (1/2), 
where the first factor of half reflects the fact that the 
number of loops is reduced by one, while the second fac- 
tor of one-half comes from the reconnection amplitude 
in Eqn [T] These rules are tabulated in Fig [2J and the 
important thing to note is that the open string can be 
treated on equal footing with (closed) loops in all cases, 
allowing one to generalize the singlet sector rules directly 
to the Stot = 1/2 sector case discussed here. 



GENERALIZATION OF THE SANDVIK-EVERTZ 
ALGORITHM 
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FIG. 3: Overlap graphs contributing to the estimator of 
{(m2)2) = E,,fc,(-P>. - l/4)(Pfci - 1/4)) with corresponding 
values of Wijki — WijWki (defined in text) given within paren- 
thesis. 



{S^} eigenbasis, and a term H^^° that is offdiagonal. In 
addition, one writes \ipi/2\) in this basis as |'!/'i/2t) — 
X]{s^} C'{Sf^}-t-|{S'f }). As in the original singlet sector 

case," each term q^-lt^Sf }t({Sf '}| n'^i Kli^f » 
generated by working in a 'spacetime' loop representation 
and using a combination of 'diagonal updates' in which 
some H^^~'^ is moved to a different bond b[^^ and 'loop 
updates' whereby each space-time loop is flipped with 
probability half; this loop update allows one to switch 
between diagonal and off-diagonal pieces of a given set of 
bond operators, while simultaneously sampling all pos- 
sible spin configurations {S'f } in the state at r = and 
T = 2m. The only difference with the original singlet 
sector case is that we now have precisely one open string 
in the space-time diagram, which connects the 'free spin 
at T — 0', i.e the unpaired spin in \ipi/2 t) to the 'free 
spin at T = 2m', i.e the unpaired spin in {ipif2 t L and 
which cannot be flipped, since the unpaired spins in all our 
Stot = 1/2, S^ot = 1/2 basis states have fixed z projec- 
tion of -1-1/2. Finally, as in the singlet case, we can easily 
generalize this procedure to treat Hamiltonians that also 
contain products — P(,-Ffc' of projectors P acting on dis- 
tinct bonds b and b' of the lattice. 



With all of this in hand, it is now easy to see that 
the singlet sector algorithm of Sandvik and collabora- 
tors [l|-|3| generalizes straightforwardly to the Stot = 1/2 
sector case for Hamiltonians of the form H = ^J2b^i>, 
where each piece — of the Hamiltonian is a projec- 
tor Pfciba acting on bond b connecting spins bi and 62- 
We start with an arbitrary Stot = 1/2 S^^f. = 1/2 state 
l^i/2t)i say |^/'i/2t) = \^o,f t)- As in the original sin- 
glet sector algorithm, we wish to stochastically sample 
(■0i/2tl(— if)^'"|'0i/2t) by sampling all possible operator 
strings Or^i in the decomposition of (— iJ)'^™ into a 
sum of products, with weight for each such string being 
proportional to (V'i/2-t-l Ot™! ^k\'4'i/2-\) ■ To do this, one 
splits each Hi, into a term H^^"^ that is diagonal in the 



ESTIMATORS 

As in the singlet sector case [H-ll], physical properties 
can be calculated by taking each space-time loop diagram 
generated by the algorithm and 'cutting it at t = m' to 
obtain the overlap diagram that represents the overlap of 

nr=l -f^h? IV'l/2t) with ('01/2t ]frZm+l ' 

Consider for instance the Neel order parameter m = 
Y.i'HiSi- Clearly, {rrix) = {my) = 0. However, (m^) 
receives contributions from sites on the open string in 
the overlap diagram, since the open string, in contrast 
to (closed) loops, has only one orientation and therefore 
cannot be flipped. More formally, we may write {A'a'j f 
\Sl\Aaf t) = Wi{A'a' t \Aaj f) and note that = 
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T]i /2 if i is part of the open string in the overlap diagram 
between {A'a'j) and {Aaj) and otherwise. We thus 
find {rriz) = where the angular brackets on the right 
denote the ensemble average over the ensemble of overlap 
diagrams generated by the modified Stot = ^/'^S^^^ — 1/2 
sector algorithm outlined above. 

We now turn to (m^) = ~ 1/4))- As noted 

earlier, whenever i and j are both in the open string or 
the same (closed) loop, the corresponding weight Wij — 
1, while Wij = 1/4 when i and j do not both belong 
to the open string or the same (closed) loop. For an 
overlap diagram with closed loops of lengths la (with a — 
1,2 .. .Ni — 1) and an open string of length 1^^ =lf, the 

latter can occur in ^ ^lalf) ways where the prime 

on the sum indicates that a = /3 is disallowed, while 
the former can occur in X^^li ways. As in the singlet 



sector case, we thus obtain 



1 v-'^' 



Sail ~ 3 Sa'/3=i ^ah, whcrc the angular brackets on 
the right indicate average over the ensemble of overlap 
diagrams generated by the algorithm. This reduces to 



^4 



ll) 



(8) 



where the angular brackets on the right again denote av- 
eraging over the ensemble of overlap diagrams generated 
by the algorithm, and the important thing to note is 
that this estimator treats the open string [a — Ni) on 
the same footing as the closed loops (a = 1, 2, ... A^; — 1). 

Finally, we consider the ground state expectation value 
of the fourth power of the Neel order parameter, i.e 
{{•m?Y). To derive the estimator for this in the Stot = 
1/2, St^t ~ 1/2, we follow Sandvik and Beach and 
write ((m2)2) = (^^^. ^^^(i^,^- _ l/4)(P,, - 1/4)). As in 
Refy, we note that the estimator for this quantity differs 
from the square of the estimator for rri? only when the ac- 
tion of Pij 'interferes' with the action of P^i, i.e when the 
actual weight W,jm = {A'a'^ f \P^JPkl\Aaf t)/(-4'a'^ t 
\Aaf t) differs from the product WijWki of the indepen- 
dent weights W,j EE {A'a'j: t \Ptj\Aaf t)/{/A'a'j: t l^o/ t 
) and Wki (defined analogously to Wij). As in the singlet 
sector case, this happens only in the two cases shown in 
Fig [31 where the difference Wijki — WijWki has been tab- 
ulated. Thus, the only new calculation needed is a count 
of the number of ways in which each of the cases Fig [3] 
(a), (b), (c), (d) arise, weighted by the corresponding val- 
ues of Wijki — WijWki. It is at this step that the open 
string needs to be treated separately, since we find that 
this count for a open string in Fig [3] (c) differs from the 
analogous count for a closed loop in Fig [3] (a) by precisely 

1 ways. On the other hand, 



one: Fig [3] (a) can arise in — ways, while Fig [3] 

4/2 

both Fig|3](b) and (d) arise in precisely 21^1'g ways (with 



(c) can arise in ^Ij 



lfs = lf) for Fig|3](d). 
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FIG. 4: Variation of the ratio of Monte Carlo estimates of 
5'z(7r, tt) and energy to exact values as a function of projection 
length: for a 3 x 5 open system (top) and a 4 x 4 periodic 
system (below). 



With all this in hand, we obtain 
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(9) 



which reduces to 
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Again, the thing to note is that the presence of the 
open string only changes the estimator by an addition 
constant — | when compared to the corresponding ex- 
pression in the singlet sector case 0. 
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FIG. 5: Boundary induced spin texture decays as ~ 1/x to 
the bulk value along a direction perpendicular to the bound- 
ary: negative differences from bulk value shown along a cut 
perpendicular to the boundary for various systemsizes. Inset 
shows the < S^{x,y) > texture on the A-sublattice sites along 
a diagonal and an edge of a 81 x 81 system. 



ILLUSTRATIVE RESULTS AND OUTLOOK 



at all sites of A sublattice for a open system with sizes 
L =. As noted by Metlitski and Sachdev [l^l, the effect 
of the boundary is to decrease the sublattice magneti- 
zation near it which is then restored to its bulk value 
away from boundary i n a power-law manner ll|, |12| . As 



demonstrated by Ref |12|, this suppression goes away as 
a power law as a function of distance \f\ from the 
edge. Using our method, we can directly calculate (5'(r)) 
in an odd by odd square lattice. On general grounds, 
one expects that (— l)'"(S'(f)) will also obey this predic- 
tion of Metlitski and Sachdev, although this quantity is 
not directly related to the usual definition of the Neel 
order parameter. With this in mind, we compare our 
results with the predictions from Ref Il2l . and find ex- 
tremely good agreement, pointing to the usefulness of 
our approach. 
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By way of illustration, we show results for a x Ly 
square lattice Heisenberg antiferromagnet with L^, Ly 
odd and open boundary conditions, and for a LxL square 
lattice Heisenberg antiferromagnet with periodic bound- 
ary conditions and L even, but with one site missing. 

To benchmark the method, we first compare the results 
foi Lx — 3, Ly — 5 open boundary condition system and 
a L — 4 period boundary condition system having one 
site missing with the corresponding exact diagonaliza- 
tion results. In Fig SI we show the dependence of the 
estimators for ground state energy, S'2(7r,7r) and as a 
function of the projection power m; this performance is 
comparable to the performance of the original algorithm 
in the singlet sector, and thus our modification provides 
a viable method to study antiferromagnets forced to have 
a Stot = 1/2 ground state due to the nature of the finite 
sample. 

With this in hand, we move on to some illustrative 
physics results. For purposes of illustration, we con- 
sider a large open boundary condition system with an 
odd number of sites. In Fig [S] shows the magnetization 
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